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Abstract. We consider the two-dimensional uniformly frustrated XY model in the 
limit of small frustration, which is equivalent to an XY system, for instance a Josephson 
junction array, in a weak uniform magnetic field applied along a direction orthogonal 
to the lattice. We show that the uniform frustration (equivalently, the magnetic field) 
destabilizes the line of fixed points which characterize the critical behaviour of the XY 
model for T < Tkt, where Tkt is the Kosterlitz-Thouless transition temperature: 
the system is paramagnetic at any temperature for sufficiently small frustration. We 
predict the critical behaviour of the correlation length and of gauge- invariant magnetic 
susceptibilities as the frustration goes to zero. These predictions are fully confirmed 
by the numerical simulations. 
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1. Introduction 

The uniformly frustrated two-dimensional (2D) XY model is defined by the lattice 
Hamiltonian 

n = -J2 Re ^ u ^* y = -^cos{e x -e y + A xy ), (i) 

{xy) (xy) 

where ip x = e l9x and U xy = e tAxy . 2D arrays of coupled Josephson junctions in a magnetic 
field are interesting physical realizations of this model In this case, the sum C(P nm ) 
of the variables A xy along the links of an elementary plaquette P nm , 

C(-Pnm) = ^(rijjnj^ji+ljm) ^(n+ljm^n+ljrra+l) 

^4(n,m+l),(n+l,m+l) J 4(n,m),(n,m+1) j (2) 

is related to the flux of an external magnetic field applied along an orthogonal direction: 
C(Pnm) = o?B /$ , where a is the lattice spacing, B is the magnetic field and 
2$ = hc/e. Hamiltonian (CQ) depends on A xy through the phases U xy and thus 
the relevant physical quantity is the product of the phases around a plaquette, i.e., 
U(P) = exp[zC(P)]. If U(P) is not 1, H is frustrated. In this paper we assume U(P) 
to be independent of the chosen plaquette, i.e., that 

U(P) = e 2 ^, (3) 

with < / < 1, independent of P. Using the invariance of the Hamiltonian under the 
transformation ip x — > ip*, it is not restrictive to take / in the interval < / < 1/2. We 
will work in a finite lattice of size L 2 with periodic boundary conditions. Therefore, we 
have 

l[U{P) = l, (4) 

p 

where the product is extended over all lattice plaquettes. This implies that fL 2 must 
be an integer. 

Hamiltonian (CQ) is invariant under the local gauge transformations 

ipx -> V x ip x , U xy -> V*U xy V y , (5) 

where V x is a phase, \V X \ = 1. Physical observables must be gauge invariant. For 
such observables, the choice of the fields A xy is irrelevant: only the value of / is 
relevant. In a finite volume, this statement is strictly true only if free boundary 
conditions are taken. If one considers periodic boundary conditions, one must also 
specify the value of exp(i A xy ) along two non-trivial lattice paths that wind around 
the lattice (they are sometimes called Polyakov loops). For instance, one must also fix 
P\{m) = exp(z2 n A( njTO)i ( n+ i im )) and P 2 {m) = exp(iJ2 n Am,n),(m,n+i)) for some fixed 
value of m. If we require the absence of magnetic circulation along these non-trivial 
paths, we must have Pi(m) = Pi(m) = 1 for any m. On a finite lattice of size L 2 , 
this condition can be satisfied only if fL is an integer, a condition that will be always 
satisfied in the numerical simulations that we shall present. 
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The critical behaviour of uniformly frustrated XY models changes dramatically 
with /. For / = the model corresponds to the standard XY model, which is not 
frustrated. It shows a Kosterlitz-Thouless transition at Tkt [on a square lattice [2] 
T KT = 0.89294(8)], where the correlation length £ diverges as ln£ ~ (T — T KT )~ X I 2 
for T > Tkt', the low-temperature phase, T < Tkt, is characterized by quasi long- 
range order — correlation functions decay algebraically — associated with a line of fixed 
points. In the case of maximal frustration, i.e. for / = 1/2, the system undergoes two 
very close continuous transitions (their critical temperature is T m 0.45 on the square 
lattice), respectively in the Ising and Kosterlitz-Thouless universality classes, see, e.g., 
[HI H] and references therein. The critical behaviour for other values of / is even more 
complex, see, e.g., [SI El El 13 U3 EH H2], and [13] for experiments. There may be 
several transitions, whose nature is not clear in most of the cases. Even the structure 
of the ground state is only partially understood [HI [TSJ EE]. For f = 1/n, where n is an 
integer number, if T c is the critical temperature where the paramagnetic phase ends, T c 
decreases with increasing n; for example, [9] T c < 0.22 for / = 1/3 and [8] T c < 0.05, 0.03 
for n = 30 and 56, respectively. These studies suggest that T c vanishes JSJ [7J as T c ~ 1 / n 
when n — > oo. The critical behaviour for irrational values of / is even less clear, see, 
e.g., [H1II2]. In this case, there are some indications that the system is paramagnetic 
for any T and that a glassy transition occurs at zero temperature |12j . 

The above-mentioned works studied the critical behaviour as a function of the 
temperature T, while keeping the uniform frustration / fixed. In this paper we 
investigate a different critical limit, i.e. we consider the limit / — >• at fixed T in the 
region T < Tkt- in other words, we investigate the effect of a small uniform frustration 
on the low-temperature XY critical behaviour. We show that a uniform frustration is 
a relevant perturbation at the fixed points that occur in the XY model for T < Tkt- 
As soon as / is non-vanishing, the correlation length becomes finite and the system is 
paramagnetic. 

The critical behaviour for small values of / can be understood within the Coulomb- 
gas picture [17J. If one considers the Villain Hamiltonian corresponding to (CQ), one can 
write the partition function as 



where [TTj Z$w is the spin-wave contribution and Hqg is the Coulomb-gas Hamiltonian: 



where is an integer (vorticity) defined at the site i of the dual lattice and V(r) is 
the lattice Coulomb potential. In ([6]) the sum over n x is restricted to configurations 
satisfying the neutrality condition [17] Y^ii n i — /) = 0. For / = and T < Tkt 
this representation allows one to show that correlations functions decay algebraically. 
The two-point correlation function is the product of a spin-wave contribution, which 
decays algebraically, and of a vortex contribution. For T < T K t charged vortices are 






(7) 
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strictly bound to form dipoles and the corresponding correlation function also decays 
algebraically [18]. For / > the picture changes. For small /, in the temperature 
interval /Tkt < T < Tkt, there are unbounded particles with n = and charge — /, 
which screen the Coulomb interaction among the vortices of charge n — f ~ n, n ^ 0. 
The Debye screening length can be easily computed. Consider a vortex of charge 1, 
surrounded by particles of charge — /. Since there is one charge — / for each lattice site, 
complete screening is achieved when these charges occupy a circle of area A, such that 
Af = 1. Thus, the screening length £ should be proportional to f^ 1 ^ 2 . In this picture, 
for / — > 0, the system is equivalent to a dilute gas (the density is proportional to f 1 ^ 2 ) of 
neutral particles interacting by means of a screened Coulomb potential V sc (r). We can 
thus perform a standard virial expansion to predict that the vortex-vortex correlation 
function is proportional to V sc (r), hence decays exponentially with a rate controlled by 
the Debye screening length. This argument indicates that, for sufficiently small / and 
any T < Tkt, the system is paramagnetic with a correlation length that scales as 

e - r 1 ' 2 , (8) 

for / 0. 

Equation (jSJ) can also be predicted by simple dimensional arguments. For a given 
value of / and T, consider a real-space renormalisation-group (RG) transformation. 
Eliminate lattice sites obtaining a lattice with a link length that is twice that of the 
original lattice. In lattice units we have £' = £/2, where we use a prime for quantities 
that refer to the decimated lattice. Analogously, we obtain /' = 4/ for the frustration 
parameter. It follows ^'f' 1 ^ 2 = if 1 / 2 - This quantity is therefore constant under RG 
transformations, i.e. £/ 1//2 = c. Under the RG transformation, the Hamiltonian 
parameters also change. In particular, the transformation induces a temperature change 
T — ► T 1 . However, for small /, one is close to the XY line of fixed points and thus 
we expect T' « T. Thus, the condition £/ 1//2 = c holds at (approximately) fixed 
temperature and / — > 0. Therefore, it implies (jHJ). 

In this paper we wish to verify numerically (J8j) and study the critical behaviour of 
gauge-invariant susceptibilities (they will be defined in the next section). Note that, in a 
sense, at fixed T < Tkt, the magnetic flux / plays the role of the reduced temperature, 
with an associated correlation- length exponent v — 1/2. 

The paper is organised as follows. In Sec. [2] we define gauge-invariant correlation 
functions, the associated susceptibilities and correlation lengths, and discuss the 
expected critical behaviour. In Sec. [3] we present some Monte Carlo (MC) results that 
fully confirm the theoretical predictions. 

2. Definitions and general scaling properties 

In order to check prediction (jSJ), we consider two different gauge- invariant correlation 
functions: 

G sq (x;y) = \{^r y )\ 2 , 
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G r (x;y) = (Re^r^]^*)- (9) 

Here T x . y is a path that connects sites x and y and [/[r^.J is a product of phases 
associated with the links that belong to T x - y . More precisely, if a link (wz) belongs to 
the path, w and z have coordinates w = (wi, w 2 ) and z — (zi, z 2 ), such that z\ — W\ > 
and z 2 — w 2 > 0, we define R wz = U wz if point w occurs before point z while moving 
along the path; otherwise, we set R wz = U^ z . The phase ^[r^] is the product of all 
the phases R wz associated with the links belonging to the path. 

The definition ([9]) of Gr(x; y) depends on a family of paths Y = {T x - y }. We assume 
this family to be translationally invariant: the path T x . y is obtained by rigidly translating 
the path T^ y _ x that connects the origin to y — x. In this case, the correlation function 
Gr{x; y) is uniquely defined by specifying the paths from the origin to any point x. 

Because of the presence of the gauge field, the Hamiltonian is not translationally 
invariant, nor is it symmetric under the symmetry transformations of the lattice. 
Nonetheless, there are generalized symmetries of the Hamiltonian that also involve gauge 
transformations. For instance, if Lf is an integer, the Hamiltonian is invariant under 
the generalized translations 

r(n,m) = V ; (n+l,m)^( nim ) i („ + l im )e " , 

i'{n ! m) — V'(n i m+l)^(n,m),(n,m+l) e27rm/ - (10) 

Gauge-invariant correlation functions are invariant under these transformations. This 
implies that they do not depend on x and y separately, but only on the difference y — x. 
This invariance can be understood intuitively if one notes that gauge-invariant quantities 
should only depend on the value of the flux through a plaquette, i.e., U(P), and of the 
Polyakov correlations Pi(m) and P 2 (m). In our model U(P) is independent of P and, if 
Lf is an integer, Pi(m) and P 2 (m) do not depend on m: hence, translation invariance 
holds. 

Analogously, the Hamiltonian is invariant under generalized transformations that 
involve lattice symmetries and gauge transformations. For instance, in infinite volume 
the Hamiltonian is invariant under the generalized reflection transformations 

|n|-l 

■0(n,m) — ^{-n,m)Km \\ [^(fc,ro),(A^l,m)^(%-l,m),(-fc,m)] ) (H) 



k=0 



where 



( m—l 



K m 



n t7 (o,*),(o,*+i) f ° r > i, 



k=0 



1 for m = 0, (12) 



-m— 1 



W U(0,k+m),(0,k+m+l) f° r TTl < —1. 

k=0 

Under these symmetries G sq (x;y) transforms covariantly. If T is a lattice symmetry, 
G sq (x;y) = G sq (Tx;Ty). These relations do not hold in general for Gr(x;y) since a 
lattice symmetry also changes the path family. 
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Given Gr(x;y) and G sq (x;y), we define the corresponding susceptibilities 

Xv = ^2G r (x;y), Xsg = ^2G sq (x;y), (13) 

y y 
where the sums are extended over all lattice points y. Because of translational invariance, 
Xsq an d Xv do not depend on the point x. Of course, xr depends on the family of paths 
T = {T^.y}. Then, for any gauge-invariant correlation function G(x;y) we define on a 
finite lattice of size L 2 

F= ^ cos[g min (?/i - xi)}G(x; y) (14) 

y=(yi,y2) 

where x = (xi,X2) and g m i n = 2n/L. The correlation length is defined by 



^ 4sin 2 (g mm /2) F ' (15) 
Note that an equally good definition of F is 

F = ^2 cos [<?min(2/2 - x 2 )]G(x;y). (16) 
y=ky\m) 

For the correlation function G sq (x; y), one can show that these two definitions of F are 
equivalent, but this is not generically the case of Gr(x;y), since this quantity is not 
symmetric under lattice transformations. In the following we use definition (I14p for F. 

In the introduction we derived a prediction for the correlation length, £ ~ f^ 1 ^ 2 - 
We wish now to obtain a similar result for the susceptibilities. In order to predict their 
scaling behaviour, let us note that, for / = and T < Tkt, ( r 4 ! o'4 } x) decays algebraically, 
i.e., (^o^x) ~ x~ t,<kT \ The critical exponent rj{T) depends on T and varies between 
r/(0) = and rj{T KT ) = 1/4. For / ^ 0, it is natural to assume that 

XT ~ / d 2 XX-^ ~ C 2 -"^ ~ r l+v(T)/2^ 

Xsq ~ / d 2 x x~ 2 ^ ~ e^ v(T) ~ f- 1+r > iT) . (17) 

In particular, these equations predict xr ~ f~ 7 ^ 8 and Xs<? ~ /~ 3//4 at T = T KT . 

The check of the previous prediction for Xsq does not present conceptual difficulties. 
Instead, when considering xr, one shoud keep in mind that this quantity depends on a 
path family. Thus, there is a natural question that should be considered first. Given a 
path family T^ 1 ^ for a given value / = f\ of the frustration parameter, we must specify 
which path family must be considered for f = fi ^ fi- Only if is chosen 
appropriately, does the relation 



-l+v(T)/2 



XT(f2) \hj 

hold for /i,/2 — > 0. A naive choice would be T^ 1 ^ = T^ 2 \ As we now discuss, this 
choice is not correct: different path families should be chosen for different values of /. 

To clarify this issue, let us imagine we are working in the continuum. For each 
/, let us consider a family of paths = {TiQ}. Because of translation invariance, 
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fl f 2 =fi/4 

r 2 Y2- 

































































































72 


























71 











































































2yi 



Figure 1. One the left we report two paths connecting the origin to y\ and 2/2, 
respectively. On the right, we report the corresponding paths connecting the origin to 
2yi and 2y 2 - The figure on the left correspond to a frustration parameter / = /1, that 
on the right to / = fi = /i/4. 



we can limit ourselves to paths going from the origin to any point y. These paths can 
be parametrised in terms of a function X^(t;y) such that X^(0;y) = for all y, 
X^\l] y) = y. The path from the origin to y is given by 

x = XV\t;y) fe[0,l]. (19) 

To determine the relation between T^ 1 ^ and T^ 2 \ one should remember that x/£ should 
be kept fixed in the critical limit. Thus, we expect the path family to be invariant only 
if all lengths are expressed in terms of £. In other words, set x = x/£f, y — y/£,f, and 
rewrite (flDl as 



x 



lxW(f;ye/) tG [0,1], 
4/ 



(20) 



where is the correlation length for the system with frustration parameter /. The 
natural requirement is therefore that the right hand side be independent of /, that is 



^X^\t ] yi h ) = ^X^\t ] yi h ). 
Since we expect £/ ~ f^ 1 ^ 2 , we obtain the relation 
X (f2 \t;ry) = rX ( - fl) (t;y), r = 



1/2 



(21) 



(22) 



In Fig. [TJwe report an example corresponding to f\ = 4/ 2 - The paths from the origin 
to yi and y 2 which belong to T^ 1 ) completely fix the paths to 2yi and 2y 2 belonging to 
r^ 2 ^. Of course, on the lattice it is impossible to ensure ff22l exactly. However, note 
that the relevant scale is fixed by the correlation length and thus, violations at the level 
of the lattice spacing are irrelevant in the critical limit. 

In the following we shall consider the path families T n = {T n -o- x }, which are specified 
by a non- negative integer n. They are defined as follows (see Fig. [2]). The path r n; o ;x 
connecting the origin to the point x = (a^,^) consists of three segments: the first one 
connects the origin to (— n, 0); the second one goes from (— n, 0) to (— n, xz); the last 
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n 



Figure 2. The path connecting the origin to the point x which belongs to the path 
family T n . 

one is horizontal, from (— n, x 2 ) to point x. We indicate with Xn{f) the corresponding 
susceptibilities and with £ n (/) the corresponding correlation lengths. These families of 
paths behave simply under the transformation (122j) . If we consider the path T n -Q- X for 
f = fi, the mapping (1221) implies that, for / = f 2 , one should consider the path T rn . . rx 
between the origin and the point rx. This implies that, if we take the path family T n for 
/ = fi, we must consider T rn for f = f 2 . As a consequence, Xn and £ n scale correctly 
only if we consider the limit n — > oo, / — > at fixed nf 1 ^ 2 . Thus, we predict the scaling 
behaviours 

Xn = r 1+ * T)/2 F x (nf 1/2 ), 

en =r 1/2 w l/2 ), (23) 

where F x (x) and F^(x) are appropriate scaling functions. In the next Section, we verify 
these predictions. 

3. Numerical results 

We perform simulations for various values of / = 1/m, m integer, and T in the interval 
T < T KT , where T KT is the critical temperature of the XY model, T KT = 0.89294(8) [2]. 
We consider finite lattices of size L 2 , where L is a multiple of 1//, and periodic boundary 
conditions for the spins. Since we perform MC simulations in a gapped phase, boundary 
conditions are expected to be irrelevant in the thermodynamic limit. Cluster algorithms 
cannot be used in the presence of frustration and thus we use an overrelaxed algorithm, 
which consists in performing microcanonical and Metropolis updates. Predictions (jHJ) 
and (1T71) hold in the thermodynamic limit, i.e. for sufficiently large values of the ratio 
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Figure 3. Scaling plot for the correlation length £ n at Tkt- For each / we report 
the data satisying n < 1/(2/). 



L/£, where finite-size effects are negligible. We find numerically that size effects are 
much smaller than our statistical errors for Lf > 3. 
In the simulations we choose the gauge 

A xy = if y = x + l, (24) 

A xy = 2-nfxi if y = x + 2, 

which is consistent with (JHJ) and with Pi(m) = P2{m) = 1, as long as L is an integer 
multiple of 1/ /. With this gauge choice the computation of the susceptibilities Xn and 
of the corresponding correlation lengths £ n is quite simple. Indeed, [/[r^.^] = 1 for any 
y if the first component of x is n, i.e., if x = (n,m), m arbitrary. Thus, if we choose 
x = (n, m) in definition ([TBI , we can compute Xn without taking into account the phases 
U xy . In practice, we have determined \n by using 

X« = ^EI]< Re VW)^>- (25) 

m y 

An analogous expression holds for the correlation lengths. 

In Figs. [3] and @] we plot the correlation lengths £ n and the susceptibilities \n at 
T = Txt for several values of / and n. In this case rj(T) = 1/4 so that Xn should scale 
as /~ 7 / 8 . It is easy to show that 

Xn = Xn+l/f, Xn = Xl/f-n, (26) 

so that in (123]) one must restrict oneself to data satisfying < n < 1/(2/). The results 
reported in the figures show the scaling behaviour (|23"1) quite precisely, confirming the 
theoretical arguments. Note that the scaling function F x (x) apparently goes to zero as 
x increases. This behaviour will be confirmed below by the analysis of a non-gauge- 
invariant correlation function. 
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Figure 4. Scaling plot for the susceptibilities \n at Tr-t- For each / we report the 
data satisying n < 1/(2/). 
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Figure 5. Critical behaviour of xo vs 1// at T = 0.4. The line is the results of a fit 
to x = a f~ e , which gives e — 0.961(3), corresponding to rj(T = 0.4) = 0.079(6). 



Good agreement is also found at T < Tkt- We check the behaviour of 
Xn=o (in this case, the same path family can be used for all values of /) up 
to T = 0.2. At T = 0.2,0.3,0.4,0.5,0.8, a fit of xo to af- 1+ ^' 2 gives 
n = 0.042(8), 0.050(6), 0.079(6), 0.098(7), 0.171(3). These results are in substantial 
agreement with the leading spin-wave contibution rj = T/(2ir), and the MC 
estimates [E] r] = 0.036(3), 0.052(5), 0.074(6), 0.100(8), 0.19(2). For example, in Fig.[5] 
we show the MC results for xo & t T = 0.4, together with the result of the fit. The data 
show a clear power-law behaviour in perfect agreement with ( fTTl) . 

We also investigated the critical behaviour of x sq , which is expected to scale as f~ 3 ^. 
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Figure 6. MC results for the non-gauge-invariant susceptibility \w and for the 
product f 3/8 Xw vs In 1// at T = T KT . 



For 1/f = 40, 60, 80, we obtain Xsq = 9.933(7), 13.630(23), 17.06(4), respectively. These 
results are fully consistent with the theoretical prediction. Indeed, the product / 3//4 Xs<? 
clearly converges to a constant as / — > (corrections are expected to be proportional 
to l/ln(l//), as in the XY model at T KT ): we have f 3/i x sq = 0.6245(5), 0.6322(11), 
0.6378(15) for the same values of /. 

Finally, we mention that correlation functions which are not gauge invariant show 
a different behaviour. For example, one may consider the susceptibility Xw associated 
with the two-point function (Re^^*) in the gauge (1251) : 

x,y 

At Tkt it shows a power-law behaviour Xw ~ f~ e a s well, but with a power e ~ 0.39, 
definitely different from the value 0.875 of the gauge-invariant definition. This result 
can be derived analytically. Indeed, we can rewrite 

1 L_1 

n=0 

where Xn is defined in (125|) . Using the properties (1261) of the susceptibilities Xn, (l28"l) can 
be rewritten as 

1/(2/) 

Xw « 2f Xn- (29) 

n=0 

In this range of values of n, as is clear from Fig. m we can use the scaling behaviour 
(T23|) and write 

-V(2/) 



fxf- 7 /* / dnF{nf^) 
Jo 
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/ .l/(2/ 1 /2) poo 

~r 3/8 / dxF(x)~f-y 8 dxF(x). (30) 

Jo Jo 

Thus, provided that F(x) is integrable (we already noted that the MC data for Xn 

are consistent with F(x) — > as x — * oo), we predict Xw ~ f~ 3 ^ 8 — f~ ' 375 , which is 

consistent with the MC data (see Fig. [6]). 

Note that the critical behaviour of Xw depends on the chosen gauge. If we use the 

gauge 

A xy = -7ifx 2 if y = x + 1, (31) 
A xy = irfxi if y = x + 2, 

the susceptibility Xw does not diverge and approaches a constant as / — ► 0. 

In conclusion, we have shown that a small amount of uniform frustration 
(equivalently, a small uniform magnetic field) destabilizes the line of fixed points that 
occur in the XY model for T < Tkt- As soon as / is different from zero, the system 
becomes paramagnetic. The critical behaviour £ ~ f^ 1 ^ 2 can be predicted by simple 
Coulomb-gas and scaling arguments. Our numerical simulations fully confirm this 
prediction. Also the scaling behaviour ffT7|) for the magnetic susceptibilities is fully 
consistent with the numerical results. 
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